more off;
pp=1;
incrlat = 20;
errlimit = 0.05;
min_num = 20;
lt_min = 12;
lt_max = 14;
kp_thrashold = 2;
%[winter - Nov,Dec,Jan; Spring - Feb,Mar,Apr;Summer - May,Jun,Jul;AUtumn -
%Aug,sep,oct];
seas_min = 1;
seas_max = 12;
load e:\manoj\projects\eej\EEJ_Stefan data;
%[iy1, im1, id1, ih1, min1, sec1] = datevec(data(:,1)+730486); % t is MJD2000;
load e:\manoj\temp\eej_season im1;
bin_data = ones([2248-216,360/incrlat])*NaN;
for i = -180:incrlat:180,
L = data(:,2) > i & data(:,2) <i+incrlat & data(:,12) > lt_min & data(:,12) < lt_max & data(:,10) <= kp_thrashold...
   & im1 >= seas_min & im1 <= seas_max;%)|im1==1);
newdata = data(L,:);
kk = find(diff(floor(newdata(:,1)))==0);
if length(kk) > 0,
    fprintf('Repeat%d\n',i);
for jk = 1:length(kk),
newdata(kk(jk)+1,6) = mean(newdata(kk(jk):kk(jk)+1,6));
end;
newdata(kk,:) = [];
clear kk;
end;
fprintf('Ndata oints %d\n',length(newdata));
kp = floor(newdata(:,1));
bin_data(kp-216,pp) = newdata(:,6);
pp=pp+1;
clear newdata;
end;

clear coefff;
clear ceff;

for i = 1:360/incrlat,
for j = 1:360/incrlat,
a = bin_data(:,i);
b = bin_data(:,j);
L = isnan(a) | isnan(b);
if sum(~L) > 2,
[c,e] = corrcoef(a(~L),b(~L));
coefff(i,j)  = c(2,1);
cerr(i,j) = e(2,1);
cnum(i,j) = sum(~L);
end;
end;
end;

X = -180+(incrlat/2):incrlat:180-(incrlat/2);
Y = -180+(incrlat/2):incrlat:180-(incrlat/2);
xx = [X(1)-(X(2)-X(1))/2 X(1:end-1)+diff(X)/2 X(end)+(X(end)-X(end-1))/2];
yy = [Y(1)-(Y(2)-Y(1))/2 Y(1:end-1)+diff(Y)/2 Y(end)+(Y(end)-Y(end-1))/2];
Data = coefff;
L = cerr < errlimit & Data < 1;
% mindata = min(min(Data(L)));
% maxdata = max(max(Data(L)));
mindata = 0;
maxdata = 1;
C = (colormap); % Reverses -> red for conductor and blue for resistor
kk=1;
for i = 1:length(X),
for j = 1:length(Y),
    if i~=j & cerr(i,j) < errlimit & cnum(i,j) > min_num,
    Xpatch(kk,:) = [xx(i) xx(i+1) xx(i+1) xx(i)];
    Ypatch(kk,:) = [yy(j) yy(j) yy(j+1) yy(j+1)];
    cidx = floor(63*(Data(i,j)-mindata)/(maxdata-mindata)+1);
    
    if cidx > 64,
        cidx = 64;
    elseif cidx < 1,
        cidx = 1;
    end;
    Cpatch(1,kk,:) = C(cidx,:);
    kk=kk+1;
    end;
end;
end;
h=patch(Xpatch',Ypatch',Cpatch);
caxis([mindata maxdata]);
colorbar;
set(h,'Edgecolor','none');
more on;